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Abstract. We carry out a numerical study of fluctuations in the 
spectrum of regular graphs. Our experiments indicate that the 
level spacing distribution of a generic fc-regular graph approaches 
Oh- that of the Gaussian Orthogonal Ensemble of random matrix the- 

| ory as we increase the number of vertices. A review of the basic 

facts on graphs and their spectra is included. 
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1. Introduction 
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O ■ A regular graph is a combinatorial structure consisting of a set V of 

O ■ 

\V\ vertices, connected by edges. Two vertices are called neighbors of 
they are connected by an edge; the graph is fc-regular if each vertex has 
exactly k neighbors. To such a graph one associates a combinatorial 
Laplacian, which operates on functions on the vertices by giving the 
sum of the differences between the values of a function / at a vertex 
and its neighbors: 



X: A/(rr) = kf(x) - V f(y) 

the sum being over all neighbors of the vertex x. The |V| eigenvalues 
= E Q < Ex < ■ ■ ■ < E\v\-i lie in the interval between and 2k. 
If we take a sequence of graphs with the number of vertices \V\ 
oo, then under certain conditions (see Section |2J) there is a limiting 
density of states analogous to Weyl's law. This gives a mean counting 
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function N(E), the expected number of levels below E, which we can 
use to measure the fluctuation properties of the eigenvalues in a large 
graph. If we "unfold" the sequence of eigenvalues (for instance by 
setting Ej = N(Ej)), then we get a sequence Ej with mean spacing 
unity: Sj := Ej + \ — Ej ~ 1. The distribution function of the spacings 
{s^ - Pjv(s) = jr Yl d( s ~ s i) ~~ is called the level spacing distribution of 
the graph. It is one of several quantities used to measure the statistical 
fluctuations of a spectrum. We wish to examine it in the limit as we 
increase the number of vertices to infinity. 

Our motivation for studying these spectral fluctuations comes from 
the theory of Quantum Chaos, where one studies fluctuations of energy 
levels of dynamical systems, for instance the spectrum of the Laplacian 
of a manifold (where the classical motion is the geodesic flow). It has 
been conjectured that generically there is a remarkable dichotomy: 

(1) If the classical dynamics are completely integrable, then 
Berry and Tabor [2| conjectured that the fluctuations are the 
same as those of an uncorrelated sequence of levels, and in par- 
ticular P(s) = e~ s is Poissonian. 

(2) If the classical dynamics are chaotic then Bohigas, Giannoni 
and Schmit |3j, jS] conjectured that the fluctuations are modeled 
by the eigenvalues of a large random symmetric matrix - the 
Gaussian Orthogonal Ensemble (GOE) 1 . 

That is, the statistics of the spectral fluctuations are universal in each 
of the two classes. 

While some obvious counter-examples exist, such as the sphere in the 
integrable case (the levels are k(k + 1) with multiplicity 2k + 1), and 
more subtle examples in the chaotic case, such as the modular surface 
(the quotient of the upper half-plane by the modular group SL(2, Z)), 
where the spacings appear to be Poissonian [I], 0, [H], |Zj, there is 
sufficient numerical evidence for us to believe that these universality 
conjectures hold in the generic case. 

1 Assuming the dynamics are invariant under time reversal. 
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In the hope of gaining some extra insight into this matter we checked 
fluctuation properties of the spectrum of a regular graph. Graphs, for 
us, will occupy an intermediate step between quantizations of genuine 
chaotic dynamical systems and the statistical models of Random Ma- 
trix Theory. While we have no direct interpretation of graphs in terms 
of classical mechanics, an analogy is the random walk on a graph: Start- 
ing with an initial probability distribution, a particle at a given vertex 
moves to one its its neighbors with equal probability. This substitute 
for dynamics is chaotic in the following sense: The walk is recurrent if 
the graph is connected (which we interpret as ergodicity), and in that 
case is mixing if the graph is not bipartite. In the bipartite case, the 
set of vertices is a union of two disjoint sets (inputs and outputs) so 
that inputs can only be connected to outputs and vice-versa. Thus if 
we start from an input vertex and walk any even number of steps then 
we will only be able to land on another input, never on an output. 

There are examples (such as some Cayley graphs, see [2], [12]) where 
there are systematic multiplicities in the spectrum and the level spac- 
ing distribution at best exists only in some singular limit. For instance 
in the case of Cayley graphs of the cyclic group Z/raZ the appropri- 
ate limit gives a rigid spectrum: E n = n, so that P(s) = 5(s — 1) 
is a Dirac delta function. Another special example, analogous to the 
modular surface, seems to have Poisson spacings (numerical evidence 
by Rockmore These examples have certain symmetries or degen- 

eracies. We tested a number of families of generic (pseudo)-random 
fc-regular graphs (see section 0] for the details of the generation algo- 
rithm). The numerical evidence we accumulated, described in Section 
indicates that the resulting family of graphs have GOE spacings. 
This should be compared with the numerical investigations by Evan- 
gelou P] which indicate that in the case of sparse random symmetric 
matrices the spacings are GOE. We are thus led to conjecture that for 
a fixed degree k > 3, the eigenvalues of the generic k- regular graph 
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on a large number of vertices have fluctuations which tend to those of 
GOE (see Section El for a more precise statement). 

The purpose of our paper is not only to describe our experimental 
results, but also to give a brief survey of the theory of Quantum Chaos 
for graph theorists, and of a bit of relevant graph theory of experts in 
Quantum Chaos. 

Accordingly, we included a survey of background material on graphs 
and their spectra in Section El, and a brief overview of Random Matrix 
Theory in Section El In section 0] we present the method used for 
generating graphs, and in section E] the results of our experiments. 

Acknowledgements. We thank N. Alon, M. Krivelevich, P. Sarnak 
and B. Sudakov for helpful conversations, and A. Odlyzko for providing 
routines to aid in the numerical computation of the GOE distribution. 
The work was partially supported by grants from the NSF, the US- 
Israel Binational Science Foundation and the Israel Science Foundation. 
D.J. was supported by an NSF postdoctoral fellowship and S.M. by an 
NSF graduate fellowhip. 

2. Graphs and their spectra 

A graph G consists of a set V of vertices and a set E of edges con- 
necting pairs of vertices. Two vertices v and w are called adjacent or 
neighboring (denoted v ~ w) if they are joined by an edge. An or- 
dering (v, w) of the endpoints of an edge e gives e an orientation; the 
second vertex is often called the head of e (denoted e+), the first one 
is called the tail (denoted e_). A graph G is directed if every edge of 
G is given an orientation. We shall mostly consider undirected graphs, 
where orientations are not specified. 

Several edges connecting the same two vertices are called multiple 
edges; a graph with multiple edges is sometimes called a multigraph 2 . 
An edge with coinciding endpoints is called a loop; we shall generally 

2 The terminology varies: occasionally what we call a graph is called a simple 
graph, while what we call a multigraph is simply called a graph. 
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consider graphs without loops or multiple edges. The degree (or va- 
lency) of a vertex is the number of edges meeting at that vertex; G is 
called k-regular if the degree of every vertex is equal to k. A walk in G 
is a sequence (vo, v i, . . . , v s ) of vertices such that Vi ~ Vi + i, it is closed 
if vq = v s . G is connected if every two vertices can be joined by a walk. 

Associated to every graph is its adjacency matrix A. It is a square 
matrix of size n = \V\ whose (i, j)-th entry is equal to the number of 
edges joining vertices Vi and Vj of G. For loopless graphs the diago- 
nal entries of A are zero. The Laplacian A is an operator acting on 
functions on the set of vertices of G. It is defined by 

(2-1) A(/)(t;) = £(/(«)-/(«;)) 

Denote by B the diagonal matrix whose i-th entry is the degree of vf, 
then 

A = B-A 

For regular graphs this gives 

(2.2) A = k -Id -A 

To motivate the analogy with the Laplace-Beltrami operator on Rie- 
mannian manifolds, we first define the incidence mapping D. To do 
that, orient all edges of G in some way. D maps functions on the set 
of vertices to functions on the set of edges by the formula 

Df(e) = f(e + ) - /(e_) 

If |V| = n and \E\ = m, the matrix of D (called the incidence matrix) 
is an n-by-m matrix whose elements are and ±1; = +1 if v j is 
the head of e,-, to —1 if it is the tail and to otherwise. The Laplacian 
matrix satisfies 

(2.3) A = DDK 

One may consider the set £ of all directed edges (|£| = 2\E\) and 
think of directed edges one of whose endpoints is v as a tangent space to 
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G at v; D can then be interpreted as a combinatorial analog of exterior 
differentiation d. The adjoint D* of D acts on functions g : £ — > R by 



Then A = £)*£), analogously to A = d*d on manifolds. 

The Laplacian is a non-negative and self-adjoint operator. A con- 
stant function on a connected component of G is an eigenfunction of 
A with eigenvalue 0; the multiplicity of is equal to the number of the 
connected components of G (exactly as for the manifold Laplacian). 
In the sequel we will only deal with connected graphs. The spectrum 
of A(G) for a Axregular graph G is clearly contained in [—k,k]; the 
spectrum of A(G) is contained in [0, 2k}. A graph is bipartite if the set 
V can be partitioned into disjoint subsets V = V\ U V 2 such that all 
edges have one endpoint in V\ and another in Vi. A /c-regular graph 
is bipartite if and only if 2k is an eigenvalue of A, and in that case 
the spectrum of A has the symmetry E i— > 2k — E. Indeed, let G be 
a bipartite graph, and let Gb be the set of the blue vertices of G, and 
G r be the set of red vertices. Let / be an eigenfunction of A(G) with 
eigenvalue E. Then let f'(v) be defined as follows: 



It is not hard to check that /' is an eigenfunction of A(G) with 
eigenvalue 2k — E. 

Denote the eigenvalues of the adjacency matrix A(G) of a A;-regular 
graph G by 



The (i, j)-th entry of the matrix A r is equal to the number of walks of 
length r starting at the vertex vi and ending at Vj. Accordingly, the 
trace of A r is equal to the number of closed walks of length r. On the 
other hand, tv(A r ) = Y^i=i K * s (by definition) equal to n times the 



D*g(v) 



E 9(e) 




v G G b 
V e G r 



k = Ai > A 2 > . . . > A n > —k 
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r-th moment of the spectral density 



(2.4) 




i=l 



of A. 



A closed walk (vq, v%, . . . , v r ) is called a cycie if v±, . . . , v r are distinct. 
The girth 7(G) of G is the length of the shortest cycle of G; all closed 
walks of length less than 7(G) necessarily involve backtracking (i.e. 
v i+i — ^i-i for some i). The number of closed walks of (necessarily 
even) length 2r < 7 starting and ending at any vertex v of a fc-regular 
graph G is equal to the number of such closed walks starting and ending 
at any vertex of the infinite /c-regular tree Tf~- 

We denote by G n ^ the set of /c-regular graphs with n vertices. It is 
known ^H] (and not hard to see) that for any fixed r > 3 the expected 
number c r (G) of r-cycles in a regular graph G G G n ^ approaches a con- 
stant as n — >• 00; accordingly, for "most" graphs G G G n ^ c r (G)/n — > 
as n — > 00. 

It is easy to show ( |14t Lemma 2.2]) that the last condition implies 
that for each fixed r and for most graphs G G G n ^ the average number 
of closed walks of length r on G is asymptotic to that of the tree. 
Accordingly, the r-th moments of the spectral density (|2.4jl approach 
those of the spectral density of the of the infinite /c-regular tree as 
n — > 00. 

It follows ([Mj) that the spectral density ()2.4j) for a general G G G n ,k 
converges to the tree density ^1] given by 



supported in Ik = \-2\J~k — 1, 2\/k — 1]. This can be regarded as an 
analog for graphs of Weyl's law for manifolds, in that both give limiting 
distributions for spectral densities. 



(2.5) 
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3. Random Matrix Theory 



We give a brief overview of the Gaussian Orthogonal Ensemble (GOE) 
of Random Matrix Theory 3 - the statistical model relevant to graphs. 
It is the space of N x N real symmetric matrices H = (Hij) with a 
probability measure P(H)dH which satisfies 

(1) P(H)dH is invariant under all orthogonal changes of basis: 

P{XHX~ 1 )dH = P(H)dH, X e O(N) 

(2) Different matrix elements are statistically independent. 
These requirements force P to be of the form 

P(H) = exp(-atr(#) 2 + btv(H) + c) 

for suitable constants a > 0, b, c. After shifting the origin and normal- 
ising one finds that the joint probability distribution of the eigenvalues 
Xj, j = 1, . . . , N of H is given by 

N 

(3.1) P N (X 1 X N )d\ = CV JJ|A; - Xj\ exp(- ^ A*) J] dX j 

i<3 j 3=1 

There is an expected limiting density for the eigenvalues of a large 
N x N matrix as iV — > oo, given by Wigner's semi-circle law: 



(3.2) R^x) 



±V2N -x 2 , \x\ < V2N 
0, \x\ > V2N 



Near the top of the semi-circle, at x — 0, the density is ^/2N fit. Thus 
if we "unfold" the eigenvalues by setting Xj := Xj \J 2N /ir, we will get 
a sequence of numbers {xj} whose mean spacing is unity, as iV — > oo. 

RMT studies spectral fluctuation of the unfolded spectrum {xj} as 
N — > oo, such the probability distribution of the nearest neighbor 
spacing s n := x n+ i —x n : For each NxN matrix H, form the probability 
measure 

1 N 

n=l 



3 The standard reference is Mehta's book |15|. 
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Then as iV — > oo, there is an expected limiting distribution P(s)ds = 
limjv^oo / p(s, H)P(H)dH called the level spacing distribution. It was 
expressed by Gaudin and Mehta in terms of a certain Fredholm deter- 
minant. For small s, P(s) ~ ^-s. 

An approximation derived by Wigner before the Gaudin-Mehta for- 
mula was known, on basis of the N = 2 case, is the Wigner surmise 

P w (s) = |«T™ V 4 

which gives a surprisingly good fit (see [TH], fig 1.5). 

It is worth emphasizing that the utility of RMT lies in that the 
predicted level spacing distribution P(s) and correlation functions are 
model-independent and appear in many instances, both probabilistic 
and deterministic, independent of features such as the level density 
(|3.2j) . For instance, numerical studies indicate that sparse random ma- 
trices have GOE spacings |9;, and the experiments described in the 
following section indicate that the same is true for eigenvalues of ran- 
dom regular graphs. 

4. Random graph generation 

We generated random fc-regular graphs using a method described 
in This method has the virtues of the ease of implementation 

and of being extremely efficient for the small (< 6) values of k of 
current interest to us. On the other hand, the running time of the 
algorithm grows exponentially with the degree k, and (at least in our 
implementation) was found impractical for k > 7 on the hardware 4 
which we used. It should be noted that in the same paper [11], Wormald 
describes an algorithm which scales well with k, but is much more 
cumbersome to implement and slower for small k. 

Wormald's algorithm is easiest explained in terms of generating ran- 
dom bipartite graphs with prescribed vertex degrees. Assume that we 
wish to generate a random bipartite graph G with blue vertices, 
named &i, . . . , &m 6 , and M r red vertices, named ri, . . . , ru r - We would 

4 A lOOMhz Pentium processor PC running Linux. 
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like the vertex h\ to have degree Vi, while the vertex rj has degree Wj. 
Evidently, we must have J^fi = J2j w j = 1-^(^)1- 

We now construct an array A of size \E(G)\. The first W\ cells of A 
contain r±, the next W2 contain rg, and so on. Now, we permute the 
E(G) cells of A by a random permutation in Se(g), to get another array 
A'. The array A' defines a bipartite (multi)graph G' as follows: The 
neighbors of b\ are the first v\ entries of A', the neighbors of 62 are the 
next t>2 entries, and so on. It is possible that G' is a multigraph, since 
two of the neighbors of some hi might well be the same. If that turns out 
to be the case, we scrap A', and generate another random permutation, 
and thus another random array A", and corresponding multigraph G", 
and so on, until we have a true bipartite graph. It is clear that if 
the valences Vi and Vj are small, this process has a good chance of 
converging in reasonable time, and it should also be intuitively fairly 
clear that each bipartite graph with prescribed degrees is equally likely 
to appear. Both statements are proved in [T9] . 

The problem of generating a random fc-regular graph can, in effect, 
be reduced to the previous problem of generating a random bipartite 
graph. To wit, to each graph G we associate a bipartite graph B G , such 
that V(B G ) = V{G) U E(G), where the blue vertices of B G correspond 
to the vertices of G, while the red vertices correspond to the edges of 
G. A vertex v is connected to e in B G , whenever e is incident to v in G. 
A fc-regular G gives rise to a graph B G , where the blue vertices have 
degree k, while the red vertices have degree 2. On the other hand, not 
every bipartite H with degrees as above arises as B G for some /c-regular 
graph G, since if H has two red vertices r± and r2 such that the blue 
neighbors of n are the same as those of r2, the corresponding G is, in 
actuality, a multigraph. 

The algorithm can thus be summarized as follows: To generate a ran- 
dom /c-regular graph with n vertices, first generate a random bipartite 
graph H with n blue vertices of degree k and nk/2 vertices of degree 2. 
If H = B G for some (obviously unique) graph G, then return G, else 
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try again. The expected running time of this process is analyzed, and 
the uniformity of the results is proved in jTHj . 

Remark. Evidently, this method is even better suited to generating 
random bipartite graphs with a prescribed degree sequence. We have 
used the algorithm to generate random 3-regular and 5-regular bipartite 
graphs. The experimental results were not substantively different from 
those for general regular graphs (as described below). 



5. Experimental results 

Once we had the adjacency matrices of the graphs constructed by the 
above method, we computed their eigenvalues. The spectral densities 
of a couple of families - one of 3-regular graphs and another of 5-regular 
graphs - are displayed in Figures 1(a) and 1(b) against McKay's law 

(Hi. 




(a) Cubic graph on 2000 vertices, (b) 5-valent graph on 500 vertices. 

Figure 1. Eigenvalue distributions of random graphs vs McKay's law 
We then unfolded the spectrum by using McKay's law, and computed 
the level spacing distribution. The resulting plots compared with GOE 
showed a good fit - see Figure 2. 




Figure 2. Level spacing distribution of a cubic graph on 2000 vertices vs GOE 

We tested the matter further by using a variant of the Kolmogorov- 
Smirnov test. One compares an empirical, sample distribution to an 
expected answer by measuring the deviation of the cumulative distri- 
bution functions of the two. Recall that if Sj, i = 1, . . . , N are random 
variables (the spacings, in our case), the empirical distribution function 
is Pjv(s) = jr ^2f = i 5(s — Si) and its cumulative distribution function is 
Cn{s) = | Si < s}. To test if the distribution function is given 

by a theoretical prediction F(s), define the discrepancy D(Cn,F) or 
Kolmogorov-Smirnov statistic to be the supremum of \Cn(s) — F(s)\ 
over s > 0. The discrepancy is small if and only if the two distributions 
are close to each other. In the case that the s« are independent, iden- 
tically distributed (definitely not the case in hand!) with cumulative 
distribution function F(s), the discrepancy goes to zero almost surely 
as N — ► oo and there is a limit law giving the the limiting distribution 
L(z) of the normalized discrepancy y/ND(C n , F) as iV — > oo: 

oo 

L(z) := lim Pr{VND(C N ,F) < z} = V (-l) j e^ 2z2 

j=—oo 

In the case that the s»'s are spacings of uncorrelated levels (hence 
certainly not independent!), the level spacing distribution is exponen- 
tial P(s) = e~ s as N —>■ oo and Pyke ^B] derives a limit law for the 
normalized discrepancy. 
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In the case where the Sj's are spacings of certain models of RMT 
(not GOE, however), Katz and Sarnak [1QJ prove that the discrepancy 
goes to zero almost surely as iV — > oo and conjecture that there is a 
limit law as in the case of Kolmogorov-Smirnov and Pyke. 

Miller (work in progress) has investigated this distribution for ran- 
dom symmetric and hermitian matrices and has numerically discovered 
that, after being normalized by multiplying by y/N, it approaches a 
limiting distribution which seems independent of the type of matrix 
involved. In Figure 3 we show this cumulative distribution function 
Lgoe( z ) of the normalized discrepancy for GOE (top plot) against 
the Kolmogorov-Smirnov "brownian bridge" L(z) (bottom plot) and 
Pyke's distribution for spacings of uncorrelated levels (middle plot). 
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Cumulative distribution functions for normalized discrepancy 

Bottom plot: Kolmogorov-Smirnov limit law. 

Middle plot: Pyke's limit law for spacings of uncorrellated levels. 

Top plot: GOE. 

Figure 3 (top and middle curves are numerical simulations). 
The numerical value of L G oe{z) can be used as a goodness-of-fit test 
to see if the eigenvalues of a large symmetric matrix have GOE spacings 
in the same way one uses the Kolmogorov-Smirnov test. 
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We computed the discrepancy for the eigenvalues of a large number of 
random graphs of particular types. Comparison of the normalized dis- 
crepancies to Miller's table gave good confidence that the spacings were 
indeed close to GOE. In Figure 4 we plot the distribution of the normal- 
ized discrepancies of a set of 4500 cubic graphs on 300 vertices against 
Miller's distribution (computed from a set of 5000 random symmetric 
120 x 120 matrices). As the figure indicates, the two distributions are 
fairly close. 



Distribution of normalized discrepancies for cubic graphs vs. GOE 



Conclusion. The numerical evidence presented above leads us to be- 
lieve that for a fixed valency k > 3, the eigenvalues of the generic 
/c-regular graph on a large number of vertices have GOE fluctuations 
in the sense that as we increase the number N of vertices, for all but 
a vanishing fraction of these graphs the discrepancy between the level 
spacing distribution of the graph and the GOE distribution goes to 
zero. 
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Figure captions 

Figure 1: Eigenvalue distribution vs McKay's law. (a) a cubic graph 
on 2000 vertices, (b) a 5-valent graph on 500 vertices. 

Figure 2. Level spacing distribution of a cubic graph on 2000 vertices 
against GOE. 

Figure 3. Cumulative distribution functions for normalized discrepency. 
Bottom plot: The Kolmogorov-Smirnov limit law. Middle plot: Pyke's 
limit law for spacings of uncorrelated levels (numerical simulation). 
Top plot: GOE (numerical simulation). 
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Figure 4. Distribution of normalized discrepancies of 4500 cubic 
graphs on 300 vertices (bar-chart) against normalized discrepancies of 
5000 random symmetric 120 x 120 matrices. 
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